In this work, we introduce a novel Quadratic Binary Optimization (QBO) framework for training a quantized neural network. The framework enables the use of arbitrary activation and loss functions through spline interpolation, while Forward Interval Propagation addresses the nonlinearities and the multi-layered, composite structure of neural networks via discretizing activation functions into linear subintervals. This preserves the universal approximation properties of neural networks while allowing complex nonlinear functions accessible to quantum solvers, broadening their applicability in artificial intelligence. Theoretically, we derive an upper bound on the approximation error and the number of Ising spins required by deriving the sample complexity of the empirical risk minimization problem from an optimization perspective. A key challenge in solving the associated large-scale Quadratic Constrained Binary Optimization (QCBO) model is the presence of numerous constraints. To overcome this, we adopt the Quantum Conditional Gradient Descent (QCGD) algorithm, which solves QCBO directly on quantum hardware. We establish the convergence of QCGD under a quantum oracle subject to randomness, bounded variance, and limited coefficient precision, and further provide an upper bound on the Time-To-Solution. To enhance scalability, we further incorporate a decomposed copositive optimization scheme that replaces the monolithic lifted model with sample-wise subproblems. This decomposition substantially reduces the quantum resource requirements and enables efficient low-bit neural network training. We further propose the usage of QCGD and Quantum Progressive Hedging (QPH) algorithm to efficiently solve the decomposed problem.

Computation with the Ising model is central to future computing technologies like quantum annealing, adiabatic quantum computing, and thermodynamic classical computing. Traditionally, computed values have been equated with ground states. This paper generalizes computation with ground states to computation with spin averages, allowing computations to take place at high temperatures. It then introduces a systematic correspondence between Ising devices and neural networks and a simple method to run trained feed-forward neural networks on Ising-type hardware. Finally, a mathematical proof is offered that these implementations are always successful.

Ising Machines are emerging hardware architectures that efficiently solve NP-Hard combinatorial optimization problems. Generally, combinatorial problems are transformed into quadratic unconstrained binary optimization (QUBO) form, but this transformation often complicates the solution landscape, degrading performance, especially for multi-state problems. To address this challenge, we model spin interactions as generalized boolean logic function to significantly reduce the exploration space. We demonstrate the effectiveness of our approach on graph coloring problem using probabilistic Ising solvers, achieving similar accuracy compared to state-of-the-art heuristics and machine learning algorithms. It also shows significant improvement over state-of-the-art QUBO-based Ising solvers, including probabilistic Ising and simulated bifurcation machines. We also design 1024-neuron all-to-all connected probabilistic Ising accelerator on FPGA with the proposed approach that shows ~10000x performance acceleration compared to GPU-based Tabucol heuristics and reducing physical neurons by 1.5-4x over baseline Ising frameworks. Thus, this work establishes superior efficiency, scalability and solution quality for multi-state optimization problems.

Photonic computing promises energy-efficient acceleration for optimization and learning, yet discrete combinatorial search and continuous function approximation have largely required distinct devices and control stacks. Here we unify k-local Ising optimization and optical Kolmogorov-Arnold network (KAN) learning on a single photonic platform, establishing a critical convergence point in optical computing. We introduce an SLM-centric primitive that realizes, in one stroke, all-optical k-local Ising interactions and fully optical KAN layers. The key idea is to convert the structural nonlinearity of a nominally linear scatterer into a per-window computational resource by adding a single relay pass through the same spatial light modulator: a folded 4f relay re-images the first Fourier plane onto the SLM so that each selected clique or channel occupies a disjoint window with its own second pass phase patch. Propagation remains linear in the optical field, yet the measured intensity in each window becomes a freely programmable polynomial of the clique sum or projection amplitude. This yields native, per clique k-local couplings without nonlinear media and, in parallel, the many independent univariate nonlinearities required by KAN layers, all trainable with in-situ physical gradients using two frames (forward and adjoint). We outline implementations on spatial photonic Ising machines, injection-locked vertical cavity surface emitting laser (VCSEL) arrays, and Microsoft analog optical computers; in all cases the hardware change is one extra lens and a fold (or an on-chip 4f loop), enabling a minimal overhead, massively parallel route to high-order Ising optimization and trainable, all-optical KAN processing on one platform.

The restricted Boltzmann machine (RBM) is a neural network based on the Ising model, well known for its ability to learn probability distributions and stochastically generate new content. However, the high computational cost of Gibbs sampling in content generation tasks imposes significant bottlenecks on electronic implementations. Here, we propose a photonic restricted Boltzmann machine (PRBM) that leverages photonic computing to accelerate Gibbs sampling, enabling efficient content generation. By introducing an efficient encoding method, the PRBM eliminates the need for computationally intensive matrix decomposition and reduces the computational complexity of Gibbs sampling from $O(N)$ to $O(1)$. Moreover, its non-Von Neumann photonic computing architecture circumvents the memory storage of interaction matrices, providing substantial advantages for large-scale RBMs. We experimentally validate the photonic-accelerated Gibbs sampling by simulating a two-dimensional Ising model, where the observed phase transition temperature closely matches the theoretical predictions. Beyond physics-inspired tasks, the PRBM demonstrates robust capabilities in generating and restoring diverse content, including images and temporal sequences, even in the presence of noise and aberrations. The scalability and reduced training cost of the PRBM framework underscore its potential as a promising pathway for advancing photonic computing in generative artificial intelligence.

Sampling Boltzmann probability distributions plays a key role in machine learning and optimization, motivating the design of hardware accelerators such as Ising machines. While the Ising model can in principle encode arbitrary optimization problems, practical implementations are often hindered by soft constraints that either slow down mixing when too strong, or fail to enforce feasibility when too weak. We introduce a two-dimensional extension of the powerful parallel tempering algorithm (PT) that addresses this challenge by adding a second dimension of replicas interpolating the penalty strengths. This scheme ensures constraint satisfaction in the final replicas, analogous to low-energy states at low temperature. The resulting two-dimensional parallel tempering algorithm (2D-PT) improves mixing in heavily constrained replicas and eliminates the need to explicitly tune the penalty strength. In a representative example of graph sparsification with copy constraints, 2D-PT achieves near-ideal mixing, with Kullback-Leibler divergence decaying as O(1/t). When applied to sparsified Wishart instances, 2D-PT yields orders of magnitude speedup over conventional PT with the same number of replicas. The method applies broadly to constrained Ising problems and can be deployed on existing Ising machines.

Oscillator Ising machines (OIMs) represent an exemplar case of using physics-inspired non-linear dynamical systems to solve computationally challenging combinatorial optimization problems (COPs). The computational performance of such systems is highly sensitive to the underlying dynamical properties, the topology of the input graph, and their relative compatibility. In this work, we explore the concept of designing different dynamical systems that minimize the same objective function but exhibit drastically different dynamical properties. Our goal is to leverage this diversification in dynamics to reduce the sensitivity of the computational performance to the underlying graph, and subsequently, enhance the overall effectiveness of such physics-based computational methods. To this end, we introduce a novel dynamical system, the Dynamical Ising Machine (DIM), which, like the OIM, minimizes the Ising Hamiltonian but offers significantly different dynamical properties. We analyze the characteristic properties of the DIM and compare them with those of the OIM. We also show that the relative performance of each model is dependent on the input graph. Our work illustrates that using multiple dynamical systems with varying properties to solve the same COP enables an effective method that is less sensitive to the input graph, while producing robust solutions.

Combinatorial optimization problems are funda- mental for various fields ranging from finance to wireless net- works. This work presents a simulated bifurcation (SB) Ising solver in CMOS for NP-hard optimization problems. Analog domain computing led to a superior implementation of this algorithm as inherent and injected noise is required in SB Ising solvers. The architecture novelties include the use of SRAM compute-in-memory (CIM) to accelerate bifurcation as well as the generation and injection of optimal decaying noise in the analog domain. We propose a novel 10-T SRAM cell capable of performing ternary multiplication. When measured with 60- node, 50% density, random, binary MAXCUT graphs, this all- to-all connected Ising solver reliably achieves above 93% of the ground state solution in 0.6us with 10.8mW average power in TSMC 180nm CMOS. Our chip achieves an order of magnitude improvement in time-to-solution and power compared to previously proposed Ising solvers in CMOS and other platforms.

Quantum or quantum-inspired Ising machines have recently shown promise in solving combinatorial optimization problems in a short time. Real-world applications, such as time division multiple access (TDMA) scheduling for wireless multi-hop networks and financial trading, require solving those problems sequentially where the size and characteristics change dynamically. However, using Ising machines involves challenges to shorten system-wide latency due to the transfer of large Ising model or the cloud access and to determine the parameters for each problem. Here we show a combinatorial optimization method using embedded Ising machines, which enables solving diverse problems at high speed without runtime parameter tuning. We customize the algorithm and circuit architecture of the simulated bifurcation-based Ising machine to compress the Ising model and accelerate computation and then built a machine learning model to estimate appropriate parameters using extensive training data. In TDMA scheduling for wireless multi-hop networks, our demonstration has shown that the sophisticated system can adapt to changes in the problem and showed that it has a speed advantage over conventional methods.

Ising machines (IM) are physics-inspired alternatives to von Neumann architectures for solving hard optimization tasks. By mapping binary variables to coupled Ising spins, IMs can naturally solve unconstrained combinatorial optimization problems such as finding maximum cuts in graphs. However, despite their importance in practical applications, constrained problems remain challenging to solve for IMs that require large quadratic energy penalties to ensure the correspondence between energy ground states and constrained optimal solutions. To relax this requirement, we propose a self-adaptive IM that iteratively shapes its energy landscape using a Lagrange relaxation of constraints and avoids prior tuning of penalties. Using a probabilistic-bit (p-bit) IM emulated in software, we benchmark our algorithm with multidimensional knapsack problems (MKP) and quadratic knapsack problems (QKP), the latter being an Ising problem with linear constraints. For QKP with 300 variables, the proposed algorithm finds better solutions than state-of-the-art IMs such as Fujitsu's Digital Annealer and requires 7,500x fewer samples. Our results show that adapting the energy landscape during the search can speed up IMs for constrained optimization.